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ABSTRACT 


This thesis compared the detection ability of the spectrogram, the 1-1/2D instan¬ 
taneous power spectrum (l-l/2Dips), the bispectrum, and outer product (dyadic) repre¬ 
sentation for digitally modulated signals corrupted by additive white Gaussian noise. 

Four detection schemes were tried on noise free BPSK, QPSK, FSK, and OOK 
signals using different transform lengths. After determining the optimum transform length, 
each test signal is corrupted by additive white Gaussian noise. Different SNR levels were 
used to determine the lowest SNR level at which the message or the modulation type 
could be extracted. 

The optimal transform length was found to be the symbol duration when processing 
BPSK, OOK, and FSK via the spectrogram, the l-l/2Dips, or the bispectrum method. 
The best transform size for QPSK was half of the symbol length. For the outer product 
(dyadic) spectral representation, the best transform size was four times larger than the 
symbol length. For all processing techniques, with the exception of the outer product 
representation, the minimum detectable SNR is about 15 dB for BPSK, FSK, and OOK 
signals and about 20 dB for QPSK signals. For the outer product spectral method, these 
values tend to be about 10 dB lower. 
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I. INTRODUCTION 


A. OVERVIEW 

The estimate of the power spectral density of discrete-time deterministic or 
stochastic signals is a useful tool in digital signal processing. In spectral 
estimation, the distribution of power as a function of frequency is obtained. The 
information contained In the power spectrum is related by the Wiener-Khintchine 
theorem to the autocorrelation sequence. This is sufficient for the complete 
statistical description of a Gaussian signal. However, one has to look beyond the 
power spectrum (autocorrelation) of a signal to extract information about 
deviations from Gaussian ity or phase relationships. 

Higher order spectra (polyspectra) defined in terms of higher order statistics 
("cumulants") of a signal, contain such Information. Particular cases of higher 
order spectra are the third and fourth-order spectra called the bispectrum and 
trispectrum, respectively. These are Fourier transforms of their respective third 
and fourth-order moments. Higher-order statistics and their spectra can be 
defined in terms of moments and cumulants. Moments and moment spectra can 
be useful in the analysis of deterministic (transient and periodic) signals whereas 
cumulants and cumulant spectra are of great importance in the analysis of 
stochastic signals [Ref. 1, 2, 9,17]. 

There are several reasons for using higher-order spectra. The first one is 
based on the property that for Gaussian signals, all cumulant spectra of order 
greater than two are identically zero. If a non-Gaussian signal is received along 
with additive Gaussian noise, the use of the higher-order cumulant will suppress 
the Gaussian noise component. 

The second reason is that polyspectra (cumulant and moment) preserve the 
true phase character of signals. This provides a technique to identify 
non-minimum phase systems or to reconstruct non-minimum phase signals. For 
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example, the Knox-Thompson and triple correlation techniques use the 
bispectrum to remove random phase distortions in images degraded by the 
propagation medium [Refs. 27, 28]. 

The third reason is based on the fact that most "real world" signals are 
non-Gaussian and thus have non-zero higher-order spectra. Hence, a signal can 
be decomposed into its higher-order spectral functions where each one of them 
contains information about the signal. This allows determination of the degree of 
deviation from Gaussianity. 

Finally, higher-order spectra are useful in analyzing the nonlinearity of a 
system operating with a random input. Polyspectra can play a key role in 
detecting and characterizing the type of nonlinearity in a system from its output. 
This allows detection and characterization of nonlinear properties in signals or 
systems. 

In this work, the spectrogram, the bispectrum, the 1-1/2D instantaneous 
power spectrum, and outer product (dyadic form) spectrum were studied with 
respect to detection and classification performance. Four detection scheme were 
tried on noise free binary phase shift keyed (BPSK), quadrature phase shift 
keyed (QPSK), frequency shift keyed (FSK), and on-off keyed (OOK) signals 
using different transform lengths. After determining the optimum transform 
length, each test signal was corrupted by additive white Gaussian noise. 
Different SNR levels were used to determine the lowest SNR level at which the 
message or the modulation type can be extracted. 

B. THESIS OUTLINE 

The thesis is presented as follows; Chapter II defines cumulants and their 
relationship to moments. Each spectral method is presented. Chapter III defines 
the generation of the test signals. Chapter IV shows how the spectrogram, the 
1-1/2D|ps, bispectrum, and the outer product representation perform on noise 
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free signals. The best transform length for each method is then determined. 
Chapter V establishes the lowest SNR level for extracting the message or the 
modulation parameters. Chapter VI presents conclusions and suggestions for 
future work. The appendices contain the Matlab programs. 
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II. PROCESSING TECHNIQUES 

This chapter describes the processing techniques used in this thesis. 
Detailed information is presented to develop an understanding of each method. 
We first define cumulants and show how they are used to form polyspectra. This 
Is followed by the spectrogram, bispectrum, 1-1/2D|ps spectrum, and the outer 
product (dyadic) form. 

A. CUMULANTS AND MOMENTS 

Given a set of n real valued random variables {xi,X 2 , .,jc„} , their joint 

moments of order r = k\ +k 2 +.... + kn are given by 


Mom[x\\x2^, .,xS"] = E{x\\x2^, . ,Xn" } 


=(-jy 


d''<^wi,w2, . ,v„) 

dw 1 * dw2 dw^....dwn" 


I 


Wi=W2=- 


.= w „=0 


( 2 . 1 ) 


where <I>(wi, W 2 , .,= ^{exp(/(>vixi + W 2 X 2 +.... + w„jc„))} is their joint 

characteristic function and E{} denotes the expectation operation. For example, 
for two random variables {x,, Xj}, we have the second-order moments Mom[x^, Xj] 
= E{Xi, X 2 }, Mom[x^^] = E{Xi^, and Mom[X 2 ^ ] = EfXj^}. 

The joint second characteristic function is defined as the natural logarithm of 
a)(wi,W 2 ,....,>v„); [Ref. 1,14] I.e., 

, W2 ,..., Wn) = ln[a>(wi, W2,...» >V„)]. (2-2) 

The joint cumulants (also called semi-invariants, [Ref.1]) of order r, of a set 
of random variables, are defined as the coefficients of the Taylor expansion of 
the second characteristic function about zero. Thus, the joint cumulants can be 
expressed in terms of the joint moments of a set of random variables. The 
cumulants are given by 


Cmw[xi’,X2%....,jc«”] =irjy 




. ,y^n) 


dwjdwn ....dwn 




I vi;i=W2=....=W»=l 




(2.3) 


5 










1. Moments and Cumulants of Stationary Processes 

If { X(k)}, k = 0,±1,±2,±3,. is a real stationary random process and its 

moments up to order n exist, then 

Mom[Xik),X{k+xx) . ,X{k+Xn-\)\ = E{X(k)X{k+x{) . X{k-¥x„-i)}. (2.4) 

will depend only on the time differences x\,x 2 .,Xn-i . [Ref. 1] We define the 

moments of a stationary random process as; 

ml{xx,X2 . ,Xn-x)^E{Xik)X(k+xx) . X{k+x„-x)} . (2.5) 

Similarly , the nth-order cumulants of {X(k)} are (n-l)-dimensional functions 
defined as 


cl(xx,X2, . ,x„.x)=Cum{X{k)X(k+xx) .JX^+t„_,)}. (2.6) 

Thus, the joint cumulants can be expressed in terms of the joint moments of 
a set of random variables. Since we only deal with a given random process, we 

simplify the notation as shown, 

mi = Mom [xj = E {xj mj = Mom [x^ ,x^] = E {x^J 

m 3 = Mom [x, ,Xi ,Xi] = E {x^} m 4 = Mom [x, ,Xi ,Xi] = E {x'’^}. 

The moments are related to the cumulants as given by 

c, = Cum [Xi] = mi 

C 2 = Cum [Xi ,xj = mj - m^^ 

C3 = Cum [Xi ,Xi ,Xi] = m3 - 3m2 m^ + 2 m®i 

C 4 = Cum [Xi ,Xi ,Xi,Xi] = m 4 - 4 m 3 m, - 3 m ^2 '• 2 m 2 m^ - 6 m^i . 

The general relationship between moments of {xi,X 2 , and joint 

cumulants Cum [x^ ,X 2 .,Xn] of order r = n is given by [Ref. 1,17] 
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Cum[x\,X2, . ,Xn\ = E(-l)^*(p - 1)!E{ n Xi}E{ IT x, n x,} , 

iGS\ i€S 2 iesp 

(2.7) 

where the summation extends over all partitions (s^.Sj.Sp), p = 1,2 .,n, of 

the set of integers (1,2,...,n). The first order cumulant is equal to the first order 
moment since only one unique partition exists for one random variable. The first 
three orders of cumulants, expressed in terms of the component indices of a set 
of random variables, are 

Ci = Cum [xj = E {xj = m^ , 

C 2 = Cum [Xi .Xj] = E{Xi, Xj} E {x,} E {Xj} = m^ 2 - m, m 2 , 

C 3 = Cum [Xi ,X 2 ,X 3 ] = E{Xi ,X 2 ,X 3 } + 2E {Xi} E {X 2 }E {X 3 } - E {x^, X 2 }E {X 3 } - 
E {Xi, X 3 }E {X 2 } - E {X 2 , X 3 }E {xJ 

= mi2,3 + 2mim2m3-mi2m3-mi3m2- m 2,3 ml . 

2. Properties of Moments and Cumulants 

Inspection of the moment-cumulant expressions show that dependencies 
among random variables are removed when cumulants are calculated. In fact, 
both the second and third order cumulants are zero if all the random variables 
are independent of each other. 

Another observation regarding the above moment-cumulant expressions 
applies in the commonly encountered situation where the means of the random 
variables are equal to zero. When this is true, the first order cumulant is zero 
since it equals the mean, the second order cumulant simply equals the variance, 
and the third order cumulant equals the third order moment. 

The relationship (2.7) implies that the computation of joint order r requires 
knowledge of all moments up to order r. 

Cumulants have three properties that make them more desirable than 
moments when it comes to higher order statistics: 
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Cumulants have three properties that make them more desirable than 
moments when it comes to higher order statistics: 

1. Each cumulant is independent of all lower order cumulants. 
Consequently, all cumulants of order greater than two are equal to zero for a 
Gaussian processs, as a Gaussian process is completely characterized by its 
first and second order moments. Higher order moments, on the other hand, can 
contain information about lower order moments. Hence, the cumulants of order 
greater than two, measure the non-Gaussian nature of a time series. [Ref. 2, 25] 

2. If a set of n random variables can be divided into more than one 
statistically independent subset, then the n®* order cumulant, unlike the n*® order 
moment, is equal to zero. [Ref. 2,18 (p. 330)] 

3. Unlike moments, the cumulants of the sum of two independent stationary 
random processes is equal to the sum of the corresponding cumulants of each 
process. [Ref. 2] 


B. CUMULANT SPECTRUM 

Suppose that the process {X(k)}, k = ....,- 2 ,- 1 , 0 , 1 , 2 ,.... is a real, stationary 

process with n*® order cumulant sequence . ,%„-i ) defined by (2.6). 

Assuming that the cumulant sequence satisfies the conditions 


I Y/ ^ 

I .... S Iwb. 

Tl=-00 T„_i=-00 

E .... S (l + |Ty|)|c^(Xi,.,T„_i)|< 00 . 

Ti=—00 Trt_i=—00 


( 2 . 8 ) 


for j = 1,2.,n-1, the n*^ order cumulant spectrum CJ( w^..w^., ) of {X(k)} 

exists, is continuous, and is defined as the (n-1) - dimensional Fourier transform 

of the n*® order cumulant sequence; [Ref.1, 22]. The conditions in (2.8) describe 
the usual conditions for a Fourier Transform to be well defined. The n*® order 
cumulant is thus defined: 
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.. ...’ X^l)exp{-X>VlTi + WiT 1 +.... + w„_lT;^l)} 

(2.9) 

where |>v,| <7c, fori=1,2.n-1 and \wi+W 2 + .+ w»-il < 7 U. 

In general, . \n„.^) is complex, i.e., has magnitude and phase. It can 

be written as 

.|C^(wi,.,>v„_i)|exp{/T^(M;i,. 

( 2 . 10 ) 

The cumulant spectrum is also periodic with period 27c , i.e., 

Q(wi,.,w„_i) = C^(wi +271,. ,Wn-\ +27t). 

C. POWER SPECTRUM 

The cumulant spectrum of order 2 is the power spectrum. It is given in (2.11) 
below. 


-foo 

C2(w)= E C2(T)exp{-7(wT)}, (2.11) 

T =—00 

for Iwl < 71, where cjW is the covariance sequence of {X(k)}. Equation (2.11) is 
also known as the Wiener-Khintchine theorem. 

As we stated earlier, the 2'’“ order cumulant is given by 
c 2 (ti) = /w 2 (ti) - (ni\y (covariance sequence) 


= m2(-Ti)-(w^)2 = cK-Ti) , (2.12) 

where m 2 {-'^\) is the autocorrelation sequence. Thus, we see that the 2"“ order 
cumulant sequence is the covariance function while the 2"'* order moment 

sequence is the autocorrelation function. From equations (2.11) and (2.12), we 
have 
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(2.13) 


Cl(w) = Cl(-w) 

cf(w) > 0 (real, nonnegative function) [Ref. 1, 9]. 


D. THE BISPECTRUM 

The bispectrum is defined as the two dimensional Fourier transform of the 
third order cumulant: 


^ + 00+00 

C3(wi,>V2) = ;-o 2 I C3(xi,T2)exp{-7(>viTi+IV 2 T 2 )} (2.14) 

{Zn) T1 =—CX) T] =-00 

where |>i;i|<7t, |w2|<7r, |wi+W2l<7C. 

Here c 3 (ti T 2 ) is the third - order cumulant sequence of {X(k)} described by: 

T2) = wK'C 1»1:2) - m] [W2(x 1) + /«2 (x 2) + ?«2(x 1 - X2)] + lim'iy 

(2.15) 

and W 3 (x 1 , X 2 ) is the third - order moment sequence. 

Analogous to the power spectrum, the bispectrum can also be defined with 

frequency domain quantities. Given N samples of stationary signal x(n), its 
Fourier transform is 


= S x(«)exp(-/>v«) (2.16) 

«=o 

Assuming that the third - order cumulant in (2.14) is computed with the 
conjugation scheme (zero mean case), i.e.; 
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^ 2 ) = E{x*{n)x{n + t i)x{n + T 2 )}, 


(2.17) 


then the equivalent frequency domain expression for the bispectrum is obtained 
through an extension of the periodogram of the form; 

C 3 (w 1, W2) = 1 )X(w 2 )X *(w 1 + W2) . (2.18) 

Equations (2.14) and (2.18) represent two different non - parametric methods 
that can be used to calculate the bispectrum. The approach in (2.14) is called the 
indirect method, while the approach in (2.18) is known as the direct method. 

In general, the bispectrum's region of support is a hexagon centered at the 
origin of the (w,, Wj) bispectral plane. Evaluation of the mean product of three 
Fourier amplitude terms shows that the bispectral plane exhibits certain 
symmetries. For a stationary, real, continuous time signal the expected value of 
the product of three Fourier components is: 


E{Xc{w{)Xc{W2)Xc{W2)} = Cc(wi, W2)5 (>vi + W 2 + W 3 ); (2.19) 

where the subscript c denotes continuous time quantities. 

Three properties of (2.19) determine the bispectral symmetry. First, since 
the frequency indexes must sum to zero, W 3 = - w, - Wj. Second, the frequency 
indexes in the expectation operation can be interchanged. Third, for the 
bispectrum of a real signal, conjugate symmetry results since C(-Wi, -Wj) = C(wi, Wj). 
Using the first property to express W 3 in terms of w, and Wj, the symmetry lines 
for a continuous time signal are [Ref. 9 & 20] 


Wl = W2. 
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2wi = -W 2 (from wi = vva), 


2 w 2 = -wi (from wi = W 3 ), 

Wi = -W2, 

W 2 = 0 (from wi = -w^) 

= 0 (from W2 = -W3). 

Important symmetry conditions follow from the properties of moments: 


cKx 1, T2) = C3(T2, Xl) = C3(-T2, Xi " X2) = ^3(12 - X1,“Xi) 


= ^3(^1 -X2,-X2) = C3(Ti,T2-Xi). 


As a consequence, knowing the third-order cumulants in any of the six 
sectors of Figure 2.1 enable us to find the entire third - order cumuiant 

sequence. These sectors include their boundaries so that, for example sector I is 
an infinite wedge bounded by the lines 12 = 0 and ti = t 2 ,ti > 0 . 

The definition of the bispectrum and the properties of the third - order 

cumulants give 


Cl(wuW2) = C^{W2,Wi) 


( 2 . 20 ) 


= Cf (-W2,-W\) = C 3 (->Vi -W 2 ,W 2 ) 


= Cl(wi,-Wi -W2) = C\{-W\ - W2,W\) 


= C^3(yV2,-W\ -W2). 
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Thus knowledge of the bispectrum in the triangular region 
W 2 >0, wi >W 2 , wi+w 2 < as shown in Figure 2.2 is sufficient for a complete 
description of the bispectrum. For real valued processes, the bispectrum has 12 

symmetry regions. 

Bispectral techniques have other applications, particularly in image 
restoration. For example, when an electromagnetic wave or sound wave 
propagates through a turbulent medium such as the ocean or the atmosphere, 
random phase perturbations scramble the phase front and induce intensity 
fluctuations from multipath interference. Images formed from these distorted 
waves lose their high spatial frequencies and become blurred. One technique, 
the Knox-Thompson and the related tripple correlation technique, utilize 
bispectral methods to combine a sequence of degraded images, while 
preserving the correlations between small frequency differences. Consequently, 
an ensemble average of the bispectra of several hundred degraded images 
allows the reconstruction of the original image, by computing the phase of the 
image from the statistical average of the phase differences. [Refs. 27, 28] 

E. THE 1-1/2 D INSTANTANEOUS POWER SPECTRUM 

The 1-1/2 D instantaneous power spectrum (1-1/2 D,ps) is a combination of 
the standard 1-1/2 D spectrum (1-1/2 Dstd) and the Instantaneous Power 
Spectrum (IPS). It is shown that in some respects the1-1/2 D,ps performs better 
than the conventional spectrogram. [Ref. 9, 16, 21] When used to process 
dynamic signals, 1-1/2 D,ps is observed to have a quicker rise and fall time than 
the spectrogram. In addition, the 1-1/2 D,ps method also appears to be a good at 
detecting low SNR stationary signals. [Ref. 16] 

The standard 1-1/2 D spectrum is a degenerate form of the bispectrum. 
When the first lag of the third order cumulant expression is set to zero, the third 
order cumulant expression becomes 
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T2) = £{x*(n)x(«)x(« + T2)} 


(2.21) 


= E{\x{n)\^x{n + X2)}> (2. 22) 

A one-dimensional Fourier transform of (2.21) produces the 1-1/2 Dstd 
spectrum. [Ref. 2] 

The instantaneous Power Spectrum is defined as the average of the 
derivatives of Page's running spectrum, [Ref. 9, 21] 

t 

= l J l‘^('f)exp(-y27i/iyT|^ ; (2.23) 

—00 

and Levin's backward running spectrum, [Ref. 9, 21] 

00 

= I J \s{x)QXip{-j2-Kfz)dx\^ . (2.24) 

t 

A discrete version of IPS is obtained by taking the Fourier transform of 
windowed correlation estimates formed using simple lag products: 

N-\ 

lPS{n, w) = ^ 2 {x(«)x* (« - A:) -h+ A:)}H'(A:)exp(-^wA:) ( 2 . 25 ) 

where w(k) Is a window function, and N is the length of the sampled data 
sequence [Ref. 9, 16, 21]. The 1-1/2D|ps is defined as 
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A^-1 

1 - \l2DiPs(n,w) = 5 Z{|x(«)|^}x(n + A:)}w(^)exp(-:/wA:) (2.26) 

*=o 


as given in [Ref. 16]. 


F. OUTER PRODUCT (DYADIC) FORM 

The outer product is the two dimensional periodogram of the sample at 
frequencies f, and fj [Ref. 23]. In particular, 

Po(f\Ji, m) = X(fi , m)X*(f 2 ,m) , (2.27) 

where 

X(fi,m)= Z (2.28) 

n-{m-l)k 

and M is the total number of segments, k is the transform length, and m is the 
segment number. 

Note that the entries of the two dimensional spectrogram have complex 
symmetry and that they are the complex conjugates of the each other. 

The outer product is averaged in a coherent and incoherent fashion. For 
the coherent averaging; all complex valued matrices are added and the absolute 
value of the resultant matrix is taken. 


P{f I )/2j rn)coH 


M 

E 

m=l 


Po(f\,f2,m) 


(2.29) 


For incoherent averaging; the absolute value of the matrices are added. 
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(2.30) 


M 

P{f\,f2,m)iNcoH= £ \P o(f\j2,ni)\ 

m=\ 

Note, the entries along coordinate (f|, fj) correspond to a conventional 
averaged (coherent or incoherent) periodogram. 
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HI. SIGNAL GENERATION 

In this chapter, we introduce the signals, BPSK, QPSK, OOK, and FSK that 
are generated to serve as test data. The computer plots of these signals are 
shown in this chapter. This chapter also covers the algorithms that are used in 
the generation of these signals. Matlab code is used for the generation of 
signals and the noise. 

A. BINARY - PHASE SHIFT KEYING (BPSK) 

Suppose the information is to be transmitted one bit at a time at some 
uniform rate R bits/s. The modulator maps the binary digit 0 into a waveform 
Si(t) and binary digit 1 into a waveform S 2 (t). In this manner, each bit from the 
channel encoder is distinctly encoded. We call this process a binary modulation. 

In binary data transmission, only one of two possible signals are sent during 
each bit interval T^. On the other hand. In an M-ary data transmission system, 
any one of M possible signals can be sent during each signaling Interval T [Ref. 
10, 11, 12]. The binary data transmission system is a special case of an M-ary 
data transmission system. 

Binary phase shift keying (BPSK) is one of the most common binary 
bandpass signaling techniques. The phase of a sinusoidal carrier is shifted by 0° 
or 180° relative to the reference sinusoid with a unipolar binary signal. It is 
equivalent to PM signaling with a digital waveform and is also equivalent to 
modulating a DSB-SC signal with a polar digital waveform. [Ref. 10, 12] 

A sinusoidal carrier wave of fixed amplitude Ac and fixed frequency fc is used 
to represent both symbols 1 and 0 . A typical BPSK is shown in Figure 3.1 and 
given by 

s{t)= AcCosilTzfct), for symbol 1 

s(t) = AcCos(^Tzfct+ 7c) for symbol 0. 
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The algorithm that is used for the generation of BPSK signal consists of the 
following steps: 

1. Creation of the message, (i.e., 1010) 

2. Creation of the sinusoidal carrier wave. 

3. Adding an appropriate phase shift (180° or 0°) to modulate the carrier. 
The transmitted signal is either in phase or out of phase with the carrier when 
signalling a binary 1 or 0, respectively. 

4. Expressing every symbol In a message as a group of bits, so that a group 
of bits represent a unit of information. For example, four bits represents a symbol 
of message. In this case, the message of step 1 will become 
1111000011110000. 

5. Adding noise (for analyzing noisy signals). 

B. FREQUENCY - SHIFT KEYING ( FSK) 

Frequency-shift keying (FSK), consists of shifting the frequency of a 
sinusoidal carrier from a mark frequency ( corresponding to a binary 1 ) to a 
space frequency ( corresponding to a binary 0 ) according to the baseband 
digital signal. 

In an FSK system, two sinusoidal waves of the same amplitude Ac but 
different frequencies f^ and fc are used to represent binary symbols 1 and 0, 
respectively. That is, we may express the binary FSK wave s(t) as 

s(t) = AcCOs(^'Kf\f), for symbol 1 

s(t) = AcCos(2n/ct), for symbol 0. 

The algorithm that is used for the generation of FSK signal consists of the 
following steps: 

1. Creation of the message, (i.e., 1010). 

2. Creation of the sinusoidal carrier wave at frequency fc. 
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3. Changing the frequency to for representing symbol 1, otherwise using f^ 
as the carrier frequency. 

4. Expressing every symbol in the message as a group of bits. So that a 
group of bits represents a unit of information. For example, four bits represent 
one symbol of the message. 

5. Adding noise (for analyzing noisy signals). 

A typical FSK signal is shown in Figure 3.2. 

C. ON - OFF KEYING ( OOK , ASK) 

On-off keying (OOK), also called amplitude shift keying (ASK), consists of 
keying (switching) a carrier sinusoid on or off with a unipolar binary signal. It is 
identical to unipolar binary modulation on a DSB-SC signal [Ref. 10,11,12]. 

In an ASK system, the binary symbol 1 is represented by transmitting a 
sinusoidal carrier wave of fixed amplitude A^, and fixed frequency fc for the bit 
duration Tt, seconds, whereas the binary symbol 0 is represented by switching off 
the carrier for seconds. In mathematical terms, we may express the binary 
ASK wave s(t) as: 

s{t) = AcCos(2ufct), for symbol 1 

= 0 , for symbol 0. 

The algorithm that is used for the generation of OOK signal consists the 
following steps; 

1. Creation of the message (i.e., 1010) 

2. Creation of the sinusoidal carrier wave. 

3. Multiplying the carrier by 1 or 0 for representing the symbol 1 or 0, 
respectively. 
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Figure 3.2: Frequency Shift Keying Generation 
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4. Expressing every symbol in the message as a group of bits, so that a 
group of bits represent a unit of information. For example, four bits represent a 
symbol of the message. In this case, the message will be 111100001 111 0000. 

5. Adding noise (for analyzing noisy signals). 

A typical OOK signal is shown in Figure 3.3. 

D. QUADRATURE PHASE (QUADRIPHASE) - SHIFT KEYING 
(QPSK) 

The quadriphase-shift keying (QPSK) is an example of M-ary data 
transmission with M = 4. In quadriphase-shift keying, one of four possible signals 
is transmitted during each signaling interval, with each signal uniquely related to 
a dibit (pairs of bits are termed dibits) [Ref. 10, 12, 20]. For example, we may 
represent the four possible dibits 00, 10, 11, and 01 (in Gray-encoded form) by 
transmitting a sinusoidal carrier with one of four possible phase values, as 
follows; 


^(r) = Acsm(2nfct + 2n) , 

for 

dibit 00 

= ^cSin(27t/c/ + 7 c/2), 

for 

dibit 10 

s(t) = AcSin(2Kfct + 7t) , 

for 

dibit 11 

s(t) = AcSm(2n/ct + 37u/2), 

for 

dibit 01 


where 0 < ^ ^ T, and Tg as symbol duration. 

Clearly, QPSK represents a special form of phase modulation. This is done 

by expressing s(t) as 

s{i) = Acsm(27zfct + (J)(r)) , 
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Figure 3.3: On - Off Keying Generation 



















where the phase (!>(/) assumes a constant value for each dibit of the incoming 
data stream. Specifically, we have 


(|)(?) = In 

for 

dibit 00, 

II 

s 

-e- 

for 

dibit 10, 

(j)(0 = 7C 

for 

dibit 11, 

KO = 371/2 

for 

dibit 01. 


The symbols can be expressed as dibits or 0, 1, 2, and 3 (or any other four 
message units). A typical QPSK signal Is shown in figure 3.4. 

The algorithm that is used for the generation of an QPSK signal consists of 
the following steps: 

1. Creation of the message that consists of four possible signals. 

2. Creation of the sinusoidal carrier wave. 

3. Adding 0° phase shift for the symbol 0, (00) (i.e., leaving it as it is), adding 
90® phase shift for the symbol 1, (10), adding 180° phase shift for the symbol 2, 
(11), adding 270° phase shift for the symbol 3, (01). 

4. Expressing every symbol in the message as a group of bits, so that a 
group of bits represent a unit of information. 

5. Adding noise (for analyzing noisy signals). 
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IV. SIMULATION 


In this chapter we compare the performance of the spectrogram, the 1-1/2 
Dips, the bispectrum, and the outer product representation. Matlab programs 
were used to conduct the computer simulations contained in this chapter. 
Simulation results are interpreted to see how each method treats the signal and 
noise components. 

For comparison BPSK, FSK, OOK, and QPSK signals were used as input 
data. Initially, simulation programs were run with signals using differently sized 
Fourier transforms. This allowed the determination of the optimum transform 
length for comparing the different methods. 


A. SIMULATION PROGRAMS 


1. SPECTRO 

The SPECTRO algorithm is given in Appendix E. It computes the 
spectrogram. The periodogram is computed as: 


Pperiy^, m) = 


{mnlen-\y\-m 

S jc(«)w(«)exp(-y>v«) 


n=m 


m = 0,\,. 


(4.1) 


where m = time index, and winlen = data length. The data segment was padded 
with zeros to the desired transform length, and w(n) was a window function. 

The user must specify the window function: 
wintype : 'O' for a rectangular window 
'T for a Hamming window. 
winlen : The desired length of data window. 

The time step was fixed at one and the spectrogram surface was calculated with 
non-normalized periodograms. 
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The time-frequency (spectrogram) surface was calculated and stored in a 
location denoted by the 'P' matrix. The columns of 'P' matrix were the frequency 
bins while the rows are the time steps. 

The spectrogram program created a series of five plots. In the upper row of 
the spectrogram figures such as Figure 4.2 the first picture was a mesh plot, the 
second one is the corresponding contour plot. The third figure is obtained by 
averaging the "P" matrix over frequency. The first picture in the lower row was 
obtained by averaging the "P" matrix over time and the second figure was the 
power spectrum of the signal using all available data. 

2. ONE_HALF D,ps 

The 1-1/2 Dips algorithm was implemented in [Ref. 16] with the extrinsic 
Matlab function ONE_HALF. Modifications were made to the original program to 
remove some undesirable smoothing steps. The revised ONE_HALF algorithm is 
contained in Appendix F. It has the same input and output format as the 
SPECTRO code. 

The user specifies an input data sequence, the size of the data window, a 
step increment, and a window function. ONE_HALF returns the time-frequency 
representation of the input sequence. The ONE_HALF algorithm was based 
upon (2.26) but was windowed to avoid undesirable product term interference 
that may form spurious peaks in the time-frequency representation. First the 
following sequence was generated: 

prod\=\element{0^element(\)—‘Slement(p^-V)\\ (4.2) 

where 
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elementQc) = |x(«)|^ {x*(« - k) +x(n + ^)}, (4.3) 

winlen was the specified data window length, and x(n) was the input data 

segment. A second sequence, designated as prod2, was formed by deleting the 
first element in prod1, adding one more element to it (i.e., element(^) ) and 
then reversing the order of the conjugated remaining elements [Ref. 91: 

prodl = [element\^),element*{^-\),....,elem^ (4.4) 

The two sequences were then concatenated with a zero inserted between them 
to form a proper correlation function having a real valued Fourier transform. The 
1-1/2 Dips estimate was obtained by multiplying the real part of the transformed 
sequence by one-half. The program computed the magnitude via the real part of 
the transform: 

1-1/2 Dips = \\ F {[prodl 0prod2]}\ (4.5) 


where F denotes the Fourier transform. 

The 1-1/2 D surface characteristics were determined by the selection of 
window type (wintype), window length (winlen) and the distance that the window 
was moved through the data sequence (step). 

The inputs are: 

data : The input observations vector. The length should be a power of 2, 

e.g., 64. 128.512. 

wintype : 'O' for Rectangular window 
'T for Hamming window 

winlen : The desired width of the window, normally half of the input length. 
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step 


: The desired time step, can be 'T or any integer (typically a power 
of 2 is used). 

The time-frequency surface was calculated and stored in a location denoted 
by the 'P' matrix. The columns of 'P' matrix were the frequency bins while the 
rows are the time steps. Time-frequency axis labeling vectors, 'freqindex' and 
'timeindex' were returned to aid in the plotting of the results. 

In the upper row of the 1-1/2D spectra figures (such as Figure 4.7) also, the 
first picture was a mesh plot, while the second one was the corresponding 
contour plot. The third figure was obtained by averaging the "P" matrix over 
frequency. The first picture in the lower row was obtained by averaging the "P" 
matrix over time and the second figure is the power spectrum of the signal using 
all available data. 

3. BISPECTRUM 

Bispectrum methods produce a frequency versus frequency representation 
instead of a time-frequency representation. The direct method for computing the 
bispectrum was realized through the Hi-Spec Matlab function BISPEC_D of [Ref. 
13]. This function returned the complex bispectrum matrix and a frequency axis 
labeling vector. The user specifies an input data sequence, the desired Fast 
Fourier Transform (FFT) size, the number of samples per segment, the amount 
of overlap between data segments, and the desired form of the smoothing 
window. The three window options available were the unit hexagonal window, 
the Parzen window, and the optimum window (also known as the minimum bias 
supremum window). 

The user must specify: 

y : data vector or time series 

nfft : FFT length 

wind : window specification for frequency-domain smoothing 
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nsamp : samples per segment 

overlap : percentage overlap 

The output of BISPEC_D.M was denoted by Bspec. The output plots of 
Bspec were an nfft x nfft array with the origin at the center, with axes pointing up 
and to the right. The another output was waxis which is a vector of frequencies 
associated with the row and column indices of Bspec. The sampling frequency 
was assumed to be 1. 

In the upper row of the bispectrum figures such as Figure 4.12, the first 
picture is the mesh plot, and the second one is the contour plot. In the lower row, 
the first picture is the frequency average plot of the ' Bspec ’ matrix , and the 
second one is the power spectrum of the signal. 

4. OUTER PRODUCT (DYADIC) FORM 

The outer product is the two dimensional periodogram of the sample at 
frequencies f^ and fz 23]. In particular, 

Po(fi Ji ,m) = X(f\,m)X*ifi.rn) , (4.6) 

where 

{km-\) 

X{fi,m)= S (4.7) 

n={m-\)k 

where M is the total number of segments, k is the transform length, and m is the 
segment number. 

Note that the entries of the two dimensional spectrogram have complex 
symmetry and that they are the complex conjugates of each other. 
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We averaged the outer product in both coherent and incoherent fashions. 
For the coherent averaging; all complex valued matrices were added and the 
absolute value of the resultant matrix was taken. 


P{f I »/2j fn)coH - 


m=\ 


(4.8) 


For incoherent averaging; the absolute value of the matrices were added. 


M 

P(f\ Jli ni)lNCOH = £ \Po(f\if2i fn)\ 


m=\ 


(4.9) 


The general procedure for averaging the outer product coherently and 
incoherently was as follows; 

1. The data, which has a length of N, was divided into segments with 
respect to the chosen transform length, so that we will have N/k segments of 
data. For example, for 256 data points and a 64 point transform length then we 
will have 4 segments of data which consist of data points 1 through 64, 65 
through 128,129 through 192, and 193 through 256. 

2. The Fourier transform of each segment was taken. 

3. Each Fourier transform was multiplied in an outer product form by its 
complex conjugate to produce k x k matrices, we have N/k matrices and each 
matrix was of size k x k. 

4. For coherent averaging; all complex valued matrices are added and the 
absolute value of the final matrix was taken (4.8). For incoherent averaging; the 
absolute value of each matrix was taken and then added (4.9). 

In the upper row of the outer product representation figures such as Figure 
4.17, the first picture is the mesh plot of the final matrix, and second one is the 
corresponding contour plot. In the lower row, the first one is a scatter plot, and 
the second one is the sum plot of the final matrix. The scatter plot is the plot of 
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all the rows (or the columns) of the final matrix. Note that the entries along 
coordinate (f|, fj) are equivalent to an averaged periodogram. In this study, the 
first quadrant of the contour plots was chosen for the analysis of contour plots 
since the contour plots had four symmetric quadrants. The outer product 
algorithm is given in Appendix G. 

B. FFT LENGTH SELECTION 

Different Fourier Transform (FFT) lengths were tried to find an optimal 
length for each of the test signal. The test message was 16 symbols long with 
each symbol of length 16, resulting in a data length of 256. In all simulations. 
Figure 4.1 through Figure 4.104 only noise free signals were used. 

1. BPSK 

The message consisted of ones and zeros. In the BPSK scheme, the 
transmitted signal was either in phase or out of phase with the carrier when 
signaling a binary 1 or 0, respectively. Generally, the sinusoidal carrier can be 
represented as follows: 

carrier s(n) = cos( 27 tj^«) 
where n = 0,1,.,255. 

In this formulation, Ni was the number of samples per symbol, and k was 
the number of cycles per symbol. Over N, samples exactly k cycles (periods) 
were produced. In this work, k was 4 and N, was 16. Hence we know the carrier 
frequency, f^, was (k / Ni)fs = (4 / 16)fs = 0.25 fg. 

This means that 16 samples represent 4 cycles. The data set was 256 
samples long, having 64 carrier cycles and representing 16 symbols. The 


35 






transform lengths which were used in this study are N, /4, Ni/2, Ni, 2Ni, and 4Ni 
where was the symbol length (i.e., = 16). 

The 16 symbol message {001 1 101 1001 1000 0} used for 
BPSK was also used for FSK, and OOK. Since every symbol was 16 data points, 
we expect a 180° phase shift at locations 33, 81, 97,129,161, 193. 

Figure 4.1 shows the plot of the sinusoidal carrier, the message, and the 
BPSK modulated signal. Figure 4.2 through Figure 4.6 shows spectrograms for 
the noise free case for FFT sizes 4, 8, 16, 32, and 64. A Hamming window was 
used. 

It can be seen that the message or the phase shifts at the relevant locations 
can be determined with window lengths of Ni /2 (which is 8), (which is 16) 
and 2 X Ni (which is 32). For the other window lengths, it was difficult to 
determine the message. For a length of 32 (in Figure 4.5), the contour plot was 
not showing, accurately, the phase shifts at their relevant locations and for a 
length of 8 (in Figure 4.3), the time averaged plot was not centered. Therefore, 
the Figure (4.4) for size 16 was the best in localizing the transition points. 

Figure 4.7 through Figure 4.11 show the 1-1/2 D spectrum for FFT sizes of 
4, 8, 16, 32, and 64. A Hamming window with a step size of 2 were used. The 
phase transition points can be seen on the mesh, contour and frequency 
averaged plots when using an FFT size of 8,16 and 32 (in Figures 4.8, 4.9, and 
4.10). For a length of 32 (in Figure 4.10), the frequency averaged plot was hard 
to interpret and for a length of 8 (in Figure 4.8), the time averaged plot was not 
centered. A size of 16 (i.e., matched signal bit length) was chosen as the 
optimum transform length for the 1-1/2 D spectral analysis since the time 
averaged plot showed the total spectral range (in Figure 4.9). 

Figure 4.12 through Figure 4.16 show the bispectrum for FFT sizes 4, 8, 16, 
32, and 64. The window type was 'T , the overlap was 50 percent, and the 
number of samples per segment (nsamp) was 4. For bispectral analysis on 
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Figure 4.1: BPSK Generation 



Figure 4.2 : BPSK Spectrogram (Window Length is 4) 








Figure 4.3: BPSK Spectrogram (Window Length is 8) 



Figure 4.4: BPSK Spectrogram (Window Length is 16) 
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Figure 4.5: BPSK Spectrogram (Window Length is 32) 



Figure 4.6: BPSK Spectrogram (Window Length is 64) 
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Figure 4.7: BPSK 1-1/2D Spectra (Window Length is 4) 



Figure 4.8: BPSK 1-1/2D Spectra (Window Length is 8) 
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Figure 4.14: BPSK Bispectrum (Fast Fourier Transform Length Is 16) 
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Figure 4.15: BPSK Bispectrum (Fast Fourier Transform Length is 32) 
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Figure 4,16; BPSK Bispectrum (Fast Fourier Transform Length is 64) 













BPSK signals, we chose an FFT length of 16 (i.e., matched signal bit length) for 
the remainder of this study. The transform length of 16 provided good frequency 
information extraction from the frequency averaged plots while adhering to the 
symbol length (see Figure 4.14). 

Figure 4.17 through Figure 4.26 are the outer product representation of 
BPSK for FFT sizes of 4, 8, 16, 32, and 64. The 4, 8, and 16 point spectra do 
not give clues nor present unique spectra about the modulation. Figures 4.23 
through 4.26 for sizes 32, and especially 64 do give clues about the modulation 
scheme. The contour plots in Figures 4.25, and 4.26 for the BPSK signal for a 
FFT size of 64 presents a picture that looks similar to the pound (#) sign which 
was very useful for differentiating BPSK from other modulation types and the 
frequency sum plots for this size presents the identical spectra for BPSK signal. 
Therefore, we choose 64 point sample for the outer product representation of the 
BPSK signal. 

2. FSK 

The simulation message consisted of ones and zeros. In the FSK scheme, 
there was no frequency change for the symbol 0, while there was a frequency 
change (or increase) for the symbol 1. 

For the FSK scheme, the modulated signal can be expressed as; 

Z(«) = sin(27c(^» , 

where a = 0 for symbol 0 

a = 1 for symbol 1 

k = 4 
N 1 = 16 
n = 0,1.255. 
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Figure 4.17: BPSK Outer Product Representation 
(Transform Length is 4, Coherent) 



Figure 4.18: BPSK Outer Product Representation 
(Transform Length is 4, Incoherent) 
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Figure 4.19: BPSK Outer Product Representation 
(Transform Length is 8, Coherent) 



Figure 4.20: BPSK Outer Product Representation 
(Transform Length is 8, Incoherent) 
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Figure 4,21: BPSK Outer Product Representation 
(Transform Length is 16, Coherent) 



Figure 4.22: BPSK Outer Product Representation 
(Transform Length is 16, Incoherent) 
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Figure 4.23: BPSK Outer Product Representation 
(Transform Length is 32, Coherent) 
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Figure 4,24: BPSK Outer Product Representation 
(Transform Length is 32, Incoherent) 
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Figure 4.25: BPSK Outer Product Representation 
(Transform Length is 64, Coherent) 


Figure 4.26: BPSK Outer Product Representation 
(Transform Length is 64, Incoherent) 
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This means that there was an increase of 1 cycle per message bit (i.e., 5 
cycles per symbol) for symbol 1 and no frequency change (i.e., 4 cycles per 
symbol) for symbol 0. Since every symbol had 16 data points, we expect 
frequency shifts of plus 1 cycle at the segments of the modulated signal 
between 33 and 81, 97 and 129, and 161 and 193. 

Figure 4.27 shows the plot of the sinusoidal carrier, the message, and the 
FSK modulated signal. Figure 4.28 through Figure 4.32 show spectrograms for 
FFT sizes 4, 8, 16, 32, and 64. A Hamming window was used. Frequency shifts 
at the proper locations were determined for FFT sizes N, /2 (which is 8), 
(which is 16) and 2 x Ni (which is 32). For the other window lengths, it was 
difficult to determine message signals. From Figure 4.30 the contour plot for a 
FFT size 16 shows well separated frequencies and the time averaged plot 
provided centered spectra. Therefore we choose 16 as the optimum transform 
length for FSK signal spectrogram analysis. 

Figure 4.33 through Figure 4.37 show the 1-1/2 D spectrum for FFT sizes 4, 
8, 16, 32, and 64. The window was a Hamming window and the step size was 2. 
On the figures, transition points (at the right times) can be seen on the mesh, 
contour and frequency averaged plots, for lengths of 8,16 and 32. But in Figure 
4.36 (for length of 32), the frequency averaged plot of 1-1/2 D spectra had more 
spikes than Figure 4.34 (for length of 8), the contour plot was not as clear as for 
the length 16. We choose 16 as the optimum transform length for 1-1/2 D 
spectral analysis of the FSK signal (see Figure 4.32). 

Figure 4.38 through Figure 4.42 show the bispectrum for FFT sizes 4, 8,16, 
32, and 64. The window type was 'T, the overlap was 50 percent, and the 
number of samples per segment (nsamp) was 16. For the bispectral analysis on 
FSK signals, an FFT length of 16 (in Figure 4.40) was chosen for the remainder 
of this study to adhere to the symbol length and to conform with the choice of the 
other spectral methods. 
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Figure 4.31: FSK Spectrogram (Window Length is 32) 
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Figure 4.32: FSK Spectrogram (Window Length is 64) 










Figure 4.33: FSK 1-1/2D Spectra (Window Length is 4) 



Figure 4.34: FSK 1-1/2D Spectra (Window Length Is 8) 






















































































Figure 4.41: FSK Bispectrum (Fast Fourier Transform Length is 32) 
















Figure 4.43 through Figure 4.52 show the outer product representation of 
the FSK signal for FFT sizes of 4, 8, 16, 32, and 64. Figures for sizes of 4, 8, and 
16 do not give clues about the modulation scheme nor present sufficiently 
different unique spectra for differentiating FSK from other modulation schemes. 
Figures 4.49 through 4.52 for sizes 32, and especially 64 give clues about the 
modulation scheme and present unique spectra for the FSK signal. In the 
contour plots and in the mesh plots of Figures 4.51 and 4.52 for FFT size of 64, 
the spectra showed the frequency shifts and the two different frequencies. 
Therefore, we choose the size of 64 for the outer product representation of the 
FSK signal (see Figures 4.51 and 4.52). 

3. OOK 

In the OOK scheme, there was no signal for the symbol 0 while there was a 
sinusoidal carrier for the symbol 1. In this study, OOK modulated signal was 
formulated as; 

x{n) = sin(27r(^^)«) 

where k = 0 for symbol 0 

k = 4 for symbol 1 

N 1 = 16 

n =0,1.255. 

We expected to see the sinusoidal carrier at the segments of the signal 
between location 33 and 81, 97 and 129, and 161 and 193. Conversely, at the 
other segments, we expected to see no carrier signal. 

Figure 4.53 shows the plot of the sinusoidal carrier, message, and the OOK 
signal. Figure 4.54 through Figure 4.58 show spectrograms for FFT sizes 4, 8, 
16, 32, and 64. A Hamming window was used. It can be seen that the changes 
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Figure 4.43: FSK Outer Product Representation 
(Transform Length is 4, Coherent) 



Figure 4.44: FSK Outer Product Representation 
(Transform Length is 4, Incoherent) 
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Figure 4.45: FSK Outer Product Representation 
(Transform Length is 8, Coherent) 
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Figure 4.46: FSK Outer Product Representation 
(Transform Length is 8, Incoherent) 
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Figure 4.47: FSK Outer Product Representation 
(Transform Length is 16, Coherent) 



Figure 4.48: FSK Outer Product Representation 
(Transform Length is 16, Incoherent) 
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Figure 4.51: FSK Outer Product Representation 
(Transform Length is 64, Coherent) 



Figure 4.52: FSK Outer Product Representation 
(Transform Length is 64, Incoherent) 






























Figure 4.55: OOK Spectrogram (Window Length Is 8) 



Figure 4.56: OOK Spectrogram (Window Length is 16) 
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at the relevant locations were apparent with the window lengths of 12 (which 
is 8), Ni (which is 16) and 2 x (which is 32). For the other window lengths, it 
was difficult to determine the message by examining the mesh plot or contour 
plot. In Figure 4.55 (i.e., size 8), the time averaged plot does not present a 
complete picture since it was not centered and in the frequency averaged plot of 
Figure 4.57 (i.e., size 32) some distortions were observed which made it hard to 
see the amplitude shifts at the relevant locations. Figure 4.56 (i.e., size 16) was 
the best since the time averaged plot was centered and it showed the amplitude 
shifts at their relevant locations in the mesh, contour, and frequency averaged 
plots. 

Figure 4.59 through Figure 4.63 show the 1-1/2 D instantaneous power 
spectrum for FFT sizes 4, 8, 16, 32, and 64. A Hamming window and a step 
size of 2 were used. Amplitude transition points at the right locations can be 
seen on the mesh, contour and frequency average plots when using FFT sizes 
of 8,16 and 32 (see Figures 4.60, 4.61, and 4.62). We choose 16 (i.e., matched 
signal bit length) as the transform length for 1-1/2 D spectral analysis of the 
OOK signal since the time averaged plot was centered and the frequency 
averaged plot showed the exact locations of the amplitude shifts (see Figure 
4.61). In Figure 4.60 for a FFT size 8, the time averaged plot was not centered 
and in Figure 4.62 we observed some distortion in the mesh, contour and 
frequency averaged plot which made it hard to locate the exact points where the 
amplitude shifts occurred. 

Figure 4.64 through Figure 4.68 show the bispectrum for FFT sizes 4, 8, 
16, 32, and 64. The window type was 'T , the overlap was 50 percent, and the 
number of samples per segment (nsamp) was 16. For the bispectral analysis, we 
choose an FFT length of 16 for the remainder of this study, adhering to the 
symbol length and conforming with the transform length choice of the 
spectrogram and 1-1/2D|ps methods for OOK. 
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Figure 4.61; OOK 1-1/2D Spectra (Window Length is 16) 



Figure 4.62: OOK 1-1/2D Spectra (Window Length is 32) 





















































Figure 4.67: OOK Bispectrum (Fast Fourier Transform Length is 32) 



Figure 4.68: OOK Bispectrum (Fast Fourier Transform Length is 64) 



























Figure 4.69: OOK Outer Product Representation 
(Transform Length is 4 , Coherent) 



Figure 4.70: OOK Outer Product Representation 
(Transform Length is 4, Incoherent) 
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Figure 4J3: OOK Outer Product Representation 
(Transform Length is 16, Coherent) 



Figure 4.74: OOK Outer Product Representation 
(Transform Length is 16, Incoherent) 
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Figure 4.77: OOK Outer Product Representation 
(Transform Length is 64, Coherent Case) 
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Figure 4.78: OOK Outer Product Representation 
(Transform Length is 64, Incoherent Case) 


















Figure 4.69 through Figure 4.78 show the outer product representation for 
FFT sizes of 4, 8, 16, 32, and 64. The figures for sizes of 4, 8, and 16 do not 
provide clues about the modulation kind. Figures 4.75 through 4.78 for sizes 32, 
and especially 64 give some information about the modulation type. We choose 
a size of 64 for the outer product representation of the OOK signal since the sum 
plot presented the classical OOK spectrum and the contour plot showed an 
identical picture for the OOK signal (four solid points) which was useful for 
differentiating the OOK signal from other signals (see Figures 4.77 and 4.78). 

4. QPSK 

In the QPSK scheme, there are four phase changes for the four possible 
symbols. In this study, the QPSK modulated signal was formulated as 

x{n) = sin(27i(^)« + 0) 

where <1) = 0 for symbol 3 
O = Till for symbol 2 
O = 7U for symbol 1 
O = 3tzI2 for symbol 0 
k = 4 

Ni = 16 

n = 0,1,.255. 

The 16 symbol message for QPSK was {3 31102103311323 
2} where every symbol was 16 data points long. We expected a 90° phase shift 
at locations 65, 97, 113, 209, 225 and 241 , a 180° phase shift at locations 33, 
81, 161, and 193 , a 270° phase shift at location 129, and no phase shifts at 
locations 17, 49, 145, and 177. The phase shifts were relative to the 
unmodulated carrier phase. 
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Figure 4.79 shows the plot of the sinusoidal carrier, message, and the 
QPSK modulated signal. Figure 4.80 through Figure 4.84 show spectrograms for 
FFT sizes of 4, 8, 16, 32, and 64. A Hamming window was used. It can be seen 
that changes at the relevant locations were revealed with window lengths of 
/2 (which is 8), N, (which is 16) and 2 x (which is 32). Figure (4.83), with 
transform length of 8, was the best since some phase transitions (at the right 
locations especially in the contour plot) were easier to see than in figures with 
other transform sizes. 

Figure 4.85 through Figure 4.89 show the 1-1/2 D spectrum for FFT sizes of 
4, 8, 16, 32, and 64. A Hamming window and the step size 2 were used. Some 
phase transition points appeared better on the mesh, contour and frequency 
averaged plots when using an FFT size of 4 and 8. We choose 8 (Figure 4.86) 
as the optimum transform length for 1-1/2 D spectra analysis of QPSK signal to 
conform with the choice made for the spectrogram. 

Figure 4.90 through Figure 4.94 show the bispectrum for FFT sizes of 4, 8, 
16, 32, and 64. The window type was 'T , overlap was 50, and number of 
samples per segment (nsamp) was 16. We choose an FFT length of 8 (Figure 
4.91) for the bispectrum analysis of the QPSK signal to conform with the choice 
made for the spectrogram and 1-1/2D|ps. 

Figure 4.95 through Figure 4.104 show the outer product representation of 
the QPSK signal for FFT sizes of 4, 8, 16, 32, and 64. The figures for FFT sizes 
of 4, 8, and 16 did not provide information about the modulation scheme. But the 
figures for sizes 32, and especially 64 showed evidence of the modulation 
scheme and present unique spectra for QPSK signal. The contour plots of 
Figures 4.103 and 4.104 for FFT size of 64 show a spot consisting of many 
peaks which represent various phase transitions and the sum plot shows the 
spike with various peaks with various heights which also represent phase shifts. 
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Figure 4.83: QPSK Spectrogram (Window Length Is 32) 
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Figure 4.84: QPSK Spectrogram {Window Length is 64) 




































Figure 4.37: QPSK 1-1/2D Spectra (Window Length is 16) 
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Figure 4.88: QPSK 1-1/2D Spectra (Window Length is 32) 













Figure 4.89; QPSK 1-1/2D Spectra (Window Length is 64) 



Figure 4.90: QPSK Bispectrum (Fast Fourier Transform Length is 4) 
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Figure 4.91: QPSK BIspectrum {Fast Fourier Transform Length is 8) 
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Figure 4.92: QPSK Bispectrum (Fast Fourier Transform Length is 16) 
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Figure 4.93; QPSK Bispectrum (Fast Fourier Transform Length is 32) 



Figure 4.94: QPSK Bispectrum (Fast Fourier Transform Length is 64) 
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Figure 4.99: QPSK Outer Product Representation 
(Transform Length is 16, Coherent) 



Figure 4.100: QPSK Outer Product Representation 
(Transform Length Is 16, Incoherent) 
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Figure 4.101: QPSK Outer Product Representation 
(Transform Length is 32, Coherent) 



Figure 4.102: QPSK Outer Product Representation 
(Transform Length is 32, Incoherent) 
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Figure 4.103: QPSK Outer Product Representation 
(Transform Length is 64, Coherent) 


Figure 4.104: QPSK Outer Product Representation 
(Transform Length is 64, Incoherent) 
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Therefore, we choose a FFT size of 64 for the outer product representation of 
the QPSK signal. 

For the remainder of this study, the following transform sizes were used: For 
BPSK, OOK, and FSK, in context of the spectrogram, the 1-1/2D,ps, and the 
bispectrum analysis, a transform length equal to the symbol length (i.e., 16) was 
chosen. Half of the symbol length (i.e., 8) was chosen for QPSK. Four times the 
symbol length (i.e., 64) was chosen as the transform length for the outer product 
representation method. 
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V. SIGNAL PLUS NOISE ANALYSIS 

In this chapter we compare the performance of the spectrogram, the 
1-1/2D,ps, the bispectrum, and the outer product representation methods. For 
the comparison, BPSK, QPSK, FSK, and OOK signals were used as input data. 
Gaussian noise was added to the original signals to create different SNR's. This 
allowed determination of the lowest SNR level at which the modulation type 
could be recognized. 

Under some conditions, the two time-frequency methods (spectrogram and 
1-1/2D|ps) could recover the original signal (message) with a sign or phase 
uncertainty. The other two frequency-frequency methods (bispectrum and outer 
product representation) only revealed the modulation type. Unless otherwise 
specified, the transform lengths for the spectrograms, the 1-1/2D|ps spectra and 
the bispectra were 16 and 64 for outer product representations. 

A. BPSK 

In this section, white Gaussian noise (AWGN) was added to the binary 
phase shift keyed (BPSK) signal. As stated, we expected 180° phase shifts at 
locations 33, 81, 97, 129, 161, and 193. For message recovery, these phase 
shifts at their appropriate locations must be detected. 

Figure 5.1 shows the plot of the sinusoidal carrier, the message, the BPSK 
signal, and the BPSK signal with noise. Figure 5.2 through Figure 5.7 show the 
spectrograms for SNR levels of infinity, 20 dB, 15 dB, 10 dB, 5 dB, and 0 dB. 

The phase shifts were apparent at the appropriate locations in the contour 
plots of all SNR levels of 15 dB or more. Transients in the contour plot represent 
phase shifts. But at the 15 dB SNR level (Figure 5.4), other transients due to 
noise appear. They corrupt the signal and confuse the observer. As a result we 
lost some ability to collect correct information about the message, modulation 
type, and bit rate for SNR levels below 20 dB. 
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Figure 5.2: BPSK Spectrogram (Infinite SNR) 
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Figure 5.3: BPSK Spectrogram ( SNR= +20 dB) 



Figure 5.4: BPSK Spectrogram ( SNR= +15 dB) 
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Figure 5.5: BPSK Spectrogram ( SNR= +10 dB) 
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Figure 5.7: BPSK Spectrogram ( SNR= +0 dB) 
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Figure 5.8 through Figure 5.13 show the 1-1/2D|ps spectra for SNR levels of 
infinity , 20 dB, 15 dB, 10 dB, 5 dB, and 0 dB. Phase shifts appear at the 
appropriate locations in the contour plots for 20 dB and above. Spreading 
(transients) in the contour plots represent these phase shifts. Below 20 dB, 
however, the noise makes the spectrum difficult to interpret. 

Figure 5.14 through Figure 5.19 show the bispectrum for SNR levels of 
infinity, 20 dB, 15 dB, 10 dB, 5 dB, and 0 dB. From the bispectrum figures, we 
were not able to see any time-related information. They do, however, represent 
unique frequency-frequency spectra which were useful to differentiate between 
modulation schemes. At the infinite SNR level, there were many peaks in the 
mesh and corresponding contour plots. The heights of the peaks were very low. 
At 20 dB and 15 dB SNR levels, the magnitudes of the peaks became higher 
with a dominant peak approximately at location (0.36, 0.36) in the contour plots. 
This peak was most likely due to noise. As the SNR level decreased to 5 dB, this 
peak was still recognizable but it was no longer the dominant peak. Also, below 
the 15 dB SNR level, some other peaks emerged. Below the 5 dB SNR level, the 
original peak at (0.36, 0.36) was not recognizable. 

Figure 5.20 through Figure 5.29 show the outer product representation 
(coherent and incoherent summation) for SNR levels of infinity, 20 dB, 10 dB, 0 
dB, and -5 dB. These figures show that the outer product form was more noise 
tolerant than the other spectral methods. This may be a result of higher 
processing gain due to the larger transform length. We were able to get 
frequency related information but not time related information. The outer product 
representation was useful for differentiating BPSK modulation from the other 
modulations by interpreting the contour and frequency sum plots. BPSK has a 
spectral signature similar to pound sign (#) in the contour plot, and two spikes 
each with almost two equal peaks in the sum plot. This information was 
recognizable until the SNR fell below 10 dB. 
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Figure 5.9: BPSK 1-1/2D Spectra (SNR= +20 dB) 



Figure 5.10: BPSK 1-1/2D Spectra (SNR= +15 dB) 
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Figure 5.16: BPSK Bispectrum (SNR= +15 dB) 


BISPECTRUM 





Figure 5.17: BPSK Bispectrum {SNR= +10 dB) 
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Figure 5.20: BPSK Outer Product Representation (Coherent, Infinite SNR) 



Figure 5.21: BPSK Outer Product Representation (Incoherent, Infinite SNR) 
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Figure 5.22: BPSK Outer Poduct Representation (Coherent, SNR= +20 dB) 



Figure 5.23: BPSK Outer Poduct Representation (Incoherent, SNR= +20 dB) 
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Figure 5.25: BPSK Outer Poduct Representation (Incoherent, SNR= +10 dB) 
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Figure 5.28: BPSK Outer Poduct Representation (Coherent, SNR= -5 dB) 


Figure 5.29: BPSK Outer Poduct Representation (Incoherent, SNR= -5 dB) 
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B. FSK 


In this section, white Gaussian noise (WGN) was added to the binary 
frequency shift keyed (FSK) signal. Different SNR levels were used to examine 
the signal In the presence of noise. As stated previously, we expected frequency 
shifts at locations 33, 81, 97, 129, 161, and 193 when interpreting 
time-frequency plots. 

Figure 5.30 shows the plot of the sinusoidal carrier, the message, the FSK 
modulated signal, and the FSK modulated signal with noise. Figure 5.31 through 
Figure 5.36 show the spectrograms for SNR levels of infinity, 20 dB, 15 dB, 10 
dB, 5 dB, and 0 dB. The frequency shifts were apparent at the appropriate 
locations in the contour plots for SNR levels of 10 dB and above. Below 10 dB, 
we were unable to recognize frequency shifts at their appropriate locations, 
despite being able to estimate, to some extent, the modulation type. As a result, 
we were loosing the ability to obtain message information, modulation type, and 
bit rate for SNR levels below 10 dB. 

Figure 5.37 through Figure 5.42 show the 1-1/2D,ps spectral results. The 
frequency shifts appeared at their appropriate locations in the contour plots for 
SNR levels of 20 dB and above. Because of the noise components, the 
spectrum became hard to interpret below the 20 dB SNR level. As a result, we 
lost the ability to get the correct information about the message, modulation type, 
and bit rate at a SNR level below 20 dB. 

Figure 5.43 through Figure 5.48 show the bispectral results. The bispectral 
figures have a unique frequency-frequency representation which was useful for 
differentiating FSK from the other modulation schemes. At an infinite SNR level, 
the FSK bispectrum pictures, especially the contour plot, were well defined (1 
high peak, and 2 smaller peaks). This information was recognizable until the 
SNR fell below 15 dB. At the low SNR levels, we saw more noise peaks (see Fig. 
5.46, 5.47 and 5.48). 
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Figure 5.30: FSK Generation With Noise 
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Figure 5.34: FSK Spectrogram (SNR= +10 dB) 
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Figure 5.38: FSK 1-1/2D Spectra (SNR= +20 dB) 



Figure 5.39: FSK 1-1/2D Spectra (SNR= +15 dB) 
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Figure 5.40: FSK 1-1/2D Spectra (SNR= +10 dB) 



Figure 5.41: FSK 1-1/2D Spectra (SNR= +5 dB) 
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Figure 5.42: FSK 1-1/2D Spectra (SNR= +0 dB) 



Figure 5.43: FSK Bispectrum (Infinite SNR) 
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Figure 5.46: FSK Bispectrum (SNR= +10 dB) 


BISPECTRUM 



Figure 5.47: FSK Bispectrum (SNR= +5 dB) 
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Figure 5.49 through Figure 5.60 show the outer product representation 
(coherent and incoherent sum) for SNR levels of infinity, 20 dB, 10 dB, 0 dB, and 
-5 dB. These figures show that the outer product form was more tolerant of noise 
than the other spectrum methods. The outer product representation was useful 
for differentiating the FSK modulation scheme from other modulation schemes 
with the aid of the contour plots and frequency sum plots. FSK outer product 
spectra had four dominant peaks in which each peak had a different height in the 
mesh plots (one high peak and three relatively shorter peaks). In the contour 
plots, these transients were symmetric across the diagonal. The first part of the 
sum plots were useful since they had two symmetric regions. FSK spectra had 
two spectral spikes of unequal height until the SNR fell below 0 dB. These 
features were unique to FSK. 

C. OOK 

In this section, white Gaussian noise (WGN) was added to the OOK signal. 
Different SNR levels were used to examine the performance of the spectral 
representations in the presence of noise. As explained in chapter III, we 
expected to see the carrier at time locations 33 to 81, 97 to 129, and 161 to 193. 
At the other locations we expected to see no carrier. For message recovery the 
amplitude shifts at the appropriate locations must be detected. 

Figure 5.61 shows the plot of the sinusoidal carrier, the message, the OOK 
modulated signal, and the OOK modulated signal with noise. Figure 5.62 through 
Figure 5.67 show the spectrograms for SNR levels of infinity, 20 dB, 15 dB, 10 
dB, 5 dB, and 0 dB. Amplitude shifts were apparent at their appropriate locations 
in the contour plots for SNR levels of 10 dB and larger. Below the 10 dB SNR 
level, we were unable to recognize the amplitude shifts, but we were still able to 
identify the modulation type. As a result we had a degraded ability to collect 
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Figure 5.49: FSK Outer Product Representation (Coherent, Infinite SNR) 



Figure 5.50: FSK Outer Product Representation (Incoherent, Infinite SNR) 




























Figure 5.53: FSK Outer Poduct Representation (Coherent, SNR= +10 dB) 



Figure 5.54: FSK Outer Poduct Representation (Incoherent, SNR= +10dB) 
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Figure 5.55: FSK Outer Poduct Representation (Coherent, SNR= +0 dB) 



Figure 5.56: FSK Outer Poduct Representation (Incoherent, SNR= +0 dB) 
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Figure 5.57: FSK Outer Poduct Representation (Coherent, SNR= -5 dB) 



Figure 5.58: FSK Outer Poduct Representation (Incoherent, SNR= -5 dB) 
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Figure 5.61: OOK Generation with Noise 



Figure 5.62: OOK Spectrogram (Infinite SNR) 
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correct information about the message, modulation type, and bit rate at and 
below 5 dB. 

Figure 5.68 through Figure 5.73 shows the 1-1/2D,ps spectra for SNR levels 
of infinity, 20 dB, 15 dB, 10 dB, 5 dB, 0 dB and -5 dB. The amplitude shifts were 
apparent at their appropriate locations for infinity, 20, 15 and 10 dB SNR levels. 
Below 10 dB, amplitude shifts did not occur at the correct locations. The 
spectrum became hard to interpret because of the noise component. As a result 
we lost the ability the message, modulation type, and bit rate for SNR levels 
below 10 dB. 

Figure 5.74 through Figure 5.79 shows the bispectrum for SNR levels of 
infinity, 20 dB, 15 dB, 10 dB, 5 dB, and 0 dB. At an infinite SNR level, the OOK 
bispectral pictures, especially the contour plots were easy to interpret. The 
magnitudes of the peaks were very low and there were many similar peaks in the 
mesh plot. The dominant peaks were located, approximately, at location 
(-0.2,-0.2; -0.2,0; 0,-0.2) in the contour plot. As the SNR level decreased, we 
saw more peaks (except at 20 dB) with larger amplitudes than the peaks at 
higher SNR levels. However, at 20 dB SNR level, the peak at location (-0.2,-0.2) 
became unrecognizable while the other peaks were still apparent and their 
magnitudes become greater with respect to the infinite SNR level. Below 20 dB 
SNR level, (until the SNR level drops to under 10 dB) the peak at location 
(-0.2,-0.2) was still identifiable but below 10 dB SNR the other two peaks were 
lost. 

Figure 5.80 through Figure 5.91 shows the outer product representation 
(coherent and incoherent sum) for SNR levels of infinity, 20 dB, 10 dB, 0 dB, and 
-5 dB. These figures show that the outer product form was more tolerant of 
noise than the other approaches. The outer product representation was useful 
for differentiating the OOK modulation scheme from other modulation schemes 
when using contour plots and frequency sum plots. OOK had one solid spike in 


136 



Figure 5.68: OOK 1-1/2D Spectra (Infinite SNR) 


1-1/2D Spectra 




Figure 5.69: OOK 1-1/2D Spectra (SNR= +20dB) 
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Figure 5.72: OOK 1-1/2D Spectra (SNR= +5dB) 



Figure 5.73: OOK 1-1/2D Spectra (SNR= 0 dB) 
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Figure 6.76: OOK Bispectrum (SNR= +15 dB) 



Figure 5.77: OOK Bispectrum (SNR= +10 dB) 
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Figure 5.80: OOK Outer Product Representation (Coherent, infinite SNR) 
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Figure 5.81: OOK Outer Product Representation (Incoherent, Infinite SNR) 




































Figure 5.82: OOK Outer Poduct Representation (Coherent, SNR= +20 dB) 



Figure 5.83: OOK Outer Poduct Representation (Incoherent, SNR= +20 dB) 
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Figure 5.86: OOK Outer Poduct Representation (Coherent, SNR= +0dB) 



Figure 5.87: OOK Outer Poduct Representation (Incoherent, SNR= +0dB) 




Figure 5.88: OOK Outer Poduct Representation (Coherent, SNR= -5dB) 



Figure 5.89: OOK Outer Poduct Representation (Incoherent, SNR= -5dB) 











Figure 5.90: OOK Outer Poduct Representation (Coherent, SNR= -10 dB) 



Figure 5.91: OOK Outer Poduct Representation (Incoherent, SNR= -10 dB) 












the contour plots. In the mesh plots, there were four dominant transients that 
was connected to each other with lines produced by smaller transients. The 
whole shape was similar to a square (this is apparent when Fig. 5.84 and Fig. 
5.86 were compared to each other). In the sum plots, OOK had one spectral 
spike. This information was distinguishable until the SNR fell below 10 dB. We 
were still able to identify the modulation type at 0 dB. These features were 
unique to OOK. 

D. QPSK 

In this section, white Gaussian noise (AWGN) was added to the quadrature 
phase shift keyed (QPSK) signal. As stated in earlier chapters, for the 
time-frequency representation we expected a 90° phase shift at locations 65, 97, 
113, 209,225 and 241 , a 180° phase shift at locations 33, 81, 161, and 193 ,and 
a 270° phase shift at location 129. The phase shifts were relative to the 
unmodulated carrier phase. The transform length for the spectrogram, the 
1-1/2D|ps spectrum, and the bispectrum was 8 while it is 64 for the outer product 
representation. 

Figure 5.92 shows a plot of the sinusoidal carrier, message, the QPSK 
modulated signal, and the QPSK modulated signal with noise. Figure 5.93 
through Figure 5.98 show the spectrograms for SNR levels of infinity, 20 dB, 15 
dB, 10 dB, 5 dB, and 0 dB. It can be seen that the phase shifts were apparent at 
their appropriate locations in contour plots for SNR levels of 20 dB and above. 
But below 20 dB, we were unable to recognize phase shift at their proper 
locations. As a result we lost the ability to identify information about the 
message, modulation type, and bit rate for SNR levels below 20 dB. 

Figure 5.99 through Figure 5.105 show the 1-1/2D,ps spectral results. Phase 
shifts were apparent at their appropriate locations in contour plots at 20 dB and 
above. Below 20 dB, these phase shifts did not occur at the correct locations. 
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Figure 5.92: QPSK Generation with Noise 
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Figure 5.93: QPSK Spectrogram (Infinite SNR) 

























Figure 5.98: QPSK Spectrogram {SNR= 0 dB) 



Figure 5.99: QPSK 1-1/2D Spectra (Infinite SNR) 
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Figure 5.103: QPSK 1-1/2D Spectra {SNR= +5dB) 


155 















Figure 5.104: QPSK 1-1/2D Spectra (SNR= OdB) 



Figure 5.105: QPSK 1-1/2D Spectra {SNR= -lOdB) 




















Because of noise components, the spectrum was difficult to interpret. Correct 
information about the message, modulation type, and bit rate for SNR levels 
below 20 dB were also difficult to determine. 

Figure 5.106 through Figure 5.111 show the bispectral results. At an infinite 
SNR level, the QPSK bispectral pictures, especially the contour and mesh plots, 
showed a peak at location approximately (0.26,0.26). As the SNR level 
decreased, the magnitude of the peaks increased. A peak was recognizable at 
20 dB and 15 dB. At SNR levels below 15 dB, this peak was not recognizable 
and any dominant peaks became harder to identify. 

Figure 5.112 through Figure 5.121 show the outer product representation 
(coherent and incoherent sum) for SNR levels of infinity, 20 dB, 10 dB, 0 dB, and 
-5 dB. From these figures, it can be seen that the outer product form was more 
tolerant to noise than the other spectrum methods. The outer product 
representation was useful for differentiating the QPSK modulation scheme from 
other modulation schemes when using the contour and frequency sum plots. In 
the contour plots, QPSK had many peaks, all together these peaks may look like 
a point. This point does not look like a solid point as in the OOK's case. In the 
sum plots, QPSK had two spectral spikes of unequal height. The width of the 
QPSK's spectral spikes is narrower than the width of FSK's spectral spikes. This 
information was identifiable until the SNR fell below 10 dB. These features were 
unique to QPSK. 
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Figure 5.108: QPSK Bispectrum (SNR= +15 dB) 



Figure 5.109: QPSK Bispectrum (SNR= +10 dB) 
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Figure 5.112: QPSK Outer Product Representation (Coherent, Infinite SNR) 



Figure 5.113: QPSK Outer Product Representation (Incoherent, Infinite SNR) 
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Figure 5.114: QPSK Outer Product Representation (Coherent, SNR= +20 dB) 



Figure 5.115: QPSK Outer Product Representation (Incoherent, SNR= +20 dB) 























Figure 5.120: QPSK Outer Product Representation (Coherent, SNR= -5 dB) 



Figure 5.121: QPSK Outer Product Representation (Incoherent, SNR= -5 dB) 
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VI. CONCLUSIONS 

A. RESULTS 

This thesis investigated four spectral analysis techniques on four signal 
modulation forms (binary phase shift keying (BPSK), frequency shift keying (FSK), 
on-off keying (OOK), and quadrature phase shift keying (QPSK)). The first part of the 
thesis determined the optimum number of points needed for each spectral analysis and 
modulation choice. The second part provided the minimum SNR level able to extract 
information about the message, modulation type, and bit rate for each spectral 
estimation method. 

1. Transform Length 

The following table summarizes the number of signal samples needed for all 
modulation types and spectral estimation methods. 



spectrogram 

1-1/2D,ps 

bispectrum 

outer product 
form 

BPSK 


Ni 

Ni 

4Ni 

QPSK 

Ni/2 

Ni/2 

Ni/2 

4Ni 

OOK 

Ni 

Ni 

Ni 

4Ni 

FSK 

Ni 

Ni 

Ni 

4Ni 


Table 6.1: Transform Length versus Modulation and Detection Type 

In this table, is the number of samples per symbol (or symbol length) which 
was 16 for all test cases. The transform length as indicated in Table 6.1 was 
determined from the ability to see time and spectral features in the spectral plots. In all 
simulations, carrier synchronization and bit synchronization were assumed. 
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2. Minimum SNR Level 


The minimum SNR levels for detecting the modulation type versus spectral 
estimation method are shown in the Table 6.2. 



spectrogram 

1-1/2D,ps 

bispectrum 

outer product 
represent. 

BPSK 

20 dB 

20 dB 

lOdB 

10 dB 

QPSK 

20 dB 

20 dB 

15dB 

10 dB 

OOK 

5dB 

10 dB 

10 dB 

OdB 

FSK 

10 dB 

20 dB 

15 dB 

OdB 


Table 6.2: The Minimum SNR Level For Identifying The Modulation Type 

Under some conditions, two methods (spectrogram and 1-1/2D|pg) could be used 
to recover the original signal (message) with some sign or phase ambiguity. The other 
two methods (bispectrum and outer product representation) could only determine the 
modulation type. 

From Table 6.2, the signal modulation techniques ranked in terms of the 
detectability, ease of classification and recognizability and arranged in decreasing order 
are: 

1. OOK 

2. FSK 

3. BPSK 

4. QPSK 

The bispectrum and outer product representations were the most noise tolerant 
methods for all modulation types investigated. Both of these techniques preserved the 
phase and frequency information of the signal. In addition the outer product 
representation had a longer transform length which improved the SNR level. The 
spectrogram was more noise tolerant than 1-1/2D|ps. The bispectrum was useful in 
identifying the modulation type. 
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B. SUGGESTIONS FOR FUTURE STUDY 

In the future, the effects of different symbol lengths, modulation types, window 
functions, and carrier and/or synchronization mismatch should be studied in the context 
of optimal processing parameters. 

Information from the bispectral plane should be examined. In particular, 
transforms along lines of constant frequency or along linear frequency indices should 
be examined for use in signal classification. 

Other cumulant spectral methods which are higher than third order could also be 
used for analyzing signals. 

Implementing an alternate bispectral calculation method relying on polar rasters 
looks like a promising starting point. Evaluation of parametric methods and comparison 
to the non-parametric methods of computing the bispectrum might prove worthwhile. 

For the outer product representation, coherent and incoherent sums should be 
studied in conjugation with different transform lengths. The phase information which 
can be obtained from the coherent sum should be analyzed. 
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APPENDIX A: BPSK CODE 


%created by Cihat YAYCI,Turkish Navy .July 1994 
%This code generates BPSK signal. 

%N is the number of information symbols (for example 16,32,64,..) 

%fs is equal to k in the formulation of sinusoidal signal cos(2*pi*k/(N*c)) 

%c is the number of bits which are used for representing each symbol. 

%d is the desired Fourier transform length, 
clear 

N=input('ENTER THE NUMBER OF INFORMATION SYMBOLS...='); 

fs =inputCENTER THE SAMPLING FREQUENCY..='): 

c=input('ENTER THE BIT LENGTH OF EACH INFORMATION SYMBOL..='): 

d=input('ENTER THE FFT LENGTH..='): 

a=ones(1:N); 

if N==16 

x16 

end 

if N==32 

x32 

end 

if N==64 

x64 

end 

ifN==128 

x128 

end 

for i=1:N; 

if x(i)<.5 

x(i)=-1: 

else 

x(i)=1: 

end 

end 

y=boxcar(c): 

data1=y*x; 

z=N*c: 

datal =reshape(data1,1 ,z) 

s=sin(2*pi*(fs/length(a)).*{1:z)): 

data=data1.*s: 

ff=fft(data): 

fff=abs(ff): 
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figure,subplot(221),plot(s) subplot{222),plot(data1) subplot(223),plot(data) 
subplot(224),plot(fl^ 
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APPENDIX B: FSK CODE 


%created by Cihat YAYCI,Turkish Navy .July 1994 
%This code generates FSK signal. 

%N is the number of information symbols (for example 16,32,64,..) 

%fs is equal to k in the formulation of sinusoidal signal cos(2*pi*k/(N*c)) 
%f is the frequency(i.e. cos(2*pi*(k+f)/(N*c)) ) 

%c is the number of bits which are used for representing each symbol 

%d is the desired Fourier transform length. 

clear N=input('ENTER THE NUMBER OF SAMPLES...='): 

fs = inputCENTER THE SAMPLING FREQUENCY..='); 

f=input(’ENTER THE FREQUENCY SHIFT...=') 

c=input(’ENTER THE DATA POINT LENGTH..='); 

d=input('ENTER THE FFT LENGTH..='): 

a=ones(1;N): 

if N==16 

x16 

end 

if N==32 

x32 

end 

if N==64 
x64 
end 

ifN==128 

x128 

end 

for i=1:N: 
if x(i)<.5 
x(i)=0; else 
x(i)=f: 
end 
end 

y=boxcar(c); 

2=N*c; 

data1=y*x: 

data 1 =reshape(data1,1 ,z) 
s=sin(2*pi*((fs+data1 )/length(a)).*(1 :z)): 
ss=fft(s); sss=abs(ss); 

figure,subplot(311),plot(data1) subplot(312),plot(s) subplot(313),plot(sss) 
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APPENDIX C: OOK CODE 

%createci by CihatYAYCI,Turkish Navy .July 1994 
%This code generates OOK signal. 

%N is the number of information symbols (for example 16,32,64,..) 

%fs is equal to k in the formulation of sinusoidal signal cos(2*pi*k/(N*c)). 
%c is the number of bits which are used for representing each symbol 
%d is the desired Fourier transform length, 
clear 

N=input('ENTER THE NUMBER OF SAMPLES...='); 

fs = inputCENTER THE SAMPLING FREQUENCY..='): 

c=input('ENTER THE DATA POINT LENGTH..=') 

d=input('ENTER THE FFT LENGTH..=’); 

a=ones(1:N): 

if N==16 

x16 

end 

if N==32 

x32 

end 

if N==64 

x64 

end 

if N==128 

x128 

end 

for i=1:N: 

if x(i)<.5 

x(i)=0: 

else 

x(i)=1: 

end 

end 

y=boxcar(c); 

data1=y*x; 

z=N*c: 

datal=reshape(data1,1 ,z) 

s=sin(2*pi*(fs/length(a)).*(1:z)): 

data=data1.*s: 

ff=fft(data); 

fff=abs(fO: 
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figure,subplot(221 ),plot(s) subplot(222),plot(data1) subplot(223),plot(data) 
subplot(224),plot(fff) 
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APPENDIX D: QPSK CODE 


%created by Cihat YAYCI,Turkish Navy .July 1994 
%This code generates QPSK signal. 

%N is the number of information symbols (for example 16,32,64,..) 

%fs is equal to k in the formulation of sinusoidal signal cos(2*pi*k/(N*c)) 
%c is the number of bits which are used for representing each symbol 
%of information symbols. 

%d is the desired Fourier transform length, 
clear 

N=input('ENTER THE NUMBER OF SAMPLES...='): 

fs = inputCENTER THE SAMPLING FREQUENCY..='): 

c=input('ENTER THE DATA POINT LENGTH..=’) 

d=input('ENTER THE FFT LENGTH..='): 

a=ones(1:N); 

if N==16 

x16 

end 

if N==32 

x32 

end 

if N==64 

x64 

end 

if N==128 

x128 

end 

fori=1:N; 

ifx(i)<.25 

x(i)=2*pi: 

end 

end 

fori=1:N; 
if x(i)<.5 
x(i)=pi/2: 
end 
end 

for i=1:N: 

ifx(i)<.75 

x(i)=pi: 

end 

end 
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fori=1:N: 
if x(i)<=1.0 
x(i)=3*pi/2: 
end 
end 

y=boxcar(c); 

data2=y*x: 

z=N*c: 

datal =reshape(data2,1 ,z) 

s1=sin((2*pi*fs/N).*(0:z-1)): 

s=sin(((2*pi*fs/N).*(0:z-1 ))+data1); 

ss=fft(s); 

sss=abs(ss): 

figure,subplot{311),plot(data1) subplot(312),plot(s) subplot(313),plot(sss) 
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APPENDIX E: SPECTROGRAM CODE 


% function [P,freqaxis,timeaxis]=spectro(data,wintype,winlen); 

% This function calculates spectrogram of the the supplied sequence, 

% "data". The user must specify: 

% "wintype" - "0" for a rectangular window. 

% "1" for a Hamming window. 

% "winlen" - The desired length of the data window. 

% The time step is fixed at one and the spectrogram is calculated with 
% non-normalized periodograms. The time-frequency surface is returned 
% in the "P" matrix. The columns of "P" are the frequency bins while the 
% rows are the time steps. Time-frequency axis labeling vectors, 

% "freqindex" and"imeindex" are also returned to aid in the plotting of results. 
% prepared by Jeff McAloon, 01 June 1993. 

% modified by Cihat Yayci, 26 July 1994. 
function [pow,freqaxis,timeaxis] = spectro(data,wintype,winlen) 
x=data; xlen=length(x): 

N=winlen; [row, col]=size(x); 
if row > 1 
x=x.': 
end 

x=[zeros(1 ,N/2) x zeros(1 ,N/2)]: 

if wintype == 1 

win=hamming(N)': 

elseif wintype == 0 

win=ones(1,N); 

end 

for ind=1:xlen 
xseg=x(ind:(N-1 +ind)): 
xseg=xseg-mean(xseg): 

SP=((abs(fft(xseg.*win,N))).''2)/1: 

POW(ind,:)=fftshift(SP): 

end 

pow=POW(:,N/2+1:N): 

[prow,pcol]=size(pow): 

%pow=POW; 

% for noise.m timeaxis=0:xlen-1; 
freqaxis=0:N/2-1; 

figure,subplot(231),mesh(pow) subplot(232),contour(pow) 
subplot(233),plot(mean(pow')) subplot(234),plot(mean(pow)) 
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APPENDIX F: ONE_HALF D,ps CODE 

%function [P,freqindex,timeindex]=one_half(data,wintype,winlen,step); 

%This function will calculate the 1 1/2 D Spectral surface. 

% The 1 1/2 D surface characteristics are determined by the selection of window 
%type (wintype), window length (winlen) and the distance that the window 
%is moved through the data sequence (step). The magnitude of the positive 
%half of the 1 1/2 D spectral plane is returned in the "P" matrix. 

%Outputs "timeindex" and "freqindex" aid in axis labeling. 

%The inputs are: 

%data -The input observations vector. The length should be a power of 2, e.g. 
%64.128,512 

%wintype: 'O' Rectangular Window 
% 'T Hamming Window 

%winlen - The desired width of the window, normally half of the input 
% length %step - Time step desired, can be 'T or a multiple of '2' 

%prepared by Karen A. Hagerman, 06 May 1992. 

%modified by Jeff McAloon, 01 June 1993. 

%modified by Cihat Yayci, 24 July 1994. 

function [P,freqindex,timeindex]=one_half(data,wintype,winlen,step) 

[datarows,datacolumns]=size(data): 

if datarows~=1 data=data.'; 

end 

siglen=length(data): 
if wintype==0 
win=ones(winlen-1,1); 
elseif 

wintype==1 

win=hamming(winlen-1): 
end 

W=[win(winlen/2;-1:1)]; 
x=[zeros(1 ,winlen) data zeros(1 ,winlen)].': 
p=zeros(siglen/step,winlen): 
index=1; 

for n=winlen+1 :step: 
siglen+winlen-step+1 
xt=x: 

il=max((winlen+1), 

(n-(winlen/2-1))); 
ul=min((n+(winlen/2-1)), 

(sigien+winlen)); 
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ctp=xt(ll:ul); 
mtp=mean(ctp); 
xt(II;ul)=ctp-mtp; 
xclt=xt((n-(winlen/2-1)): 

(n+(winlen/2-1))): 

Xn=[abs(xt(n))''2:abs(xt(n))''2]; 
Xm=[conj(xt(n:-1;n-(winlen/2-1))) xt(n:n+(winlen/2-1))]; 
product=((Xm*Xn).*W).': 

product=[product 0 conj(product(winlen/2;-1:2))]; 
pt=0.5*fftshift(fft(product)): p(index,:)=pt: 
index=index+1: 
end 

P=abs(p(:,winlen/2+1 iwinlen)); 

[Prow, Pcol]=size(P); 

%P=abs(p): 

% for noise.m 
[prow,pcolumn]=size(P); 
freqindex=[0:pcolumn-1]; timeindex=[1 :prow]: 
%subplot(221),plot(freqindex,timeindex,P) 
figure,subplot(231), 
mesh(P) 

title('1-1/2D Spectral Surface') 
xlabel('fr), ylabel('f2') 
subplot(232), 

contour{P,5,freqindex,timeindex) 

title('1-1/2D Spectral Surface') 

xlabel('fl'), ylabel('f2') 

subplot(233), 

plot(mean(P')) 

subplot(234), 

plot(mean(P)) 
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APPENDIX G: OUTER PRODUCT FORM CODE 


(for BPSK of FFT size 64) 

%created by Cihat YAYCI,Turkish Navy ,July 1994 

%This code creates the outer product representation of the BPSK signal for a FFT 
%length of 64. For the other modulated signals and different transform 
%lengths the same logic can be used. 

%N is the number of information symbols (for example 16,32,64,..) 

%fs is equal to k in the formulation of sinusoidal signal cos(2*pi*k/(N*c)) 

%f is the frequency(i.e. cos(2*pi*(k+f)/(N*c)) ) 

%c is the number of bits which are used for representing each symbol 
%d is the desired Fourier transform length. 

clear 

N=input('ENTER THE NUMBER OF SAMPLES...='); 

fs =input('ENTER THE SAMPLING FREQUENCY..='): 

c=input('ENTER THE DATA POINT LENGTH..='); 

d=input('ENTER THE FFT LENGTH..=’): 

a=ones(1:N): 

ifN==16 

x16 

end 

if N==32 

x32 

end 

if N==64 

x64 

end 

ifN==128 

x128 

end 

fori=1:N; 

ifx(i)<.5 

x(i)=-1: 

else 

x(i)=1: 

end 

end 

y=boxcar(c): 

data1=y*x: 

z=N*c; 

data 1 =reshape(data1,1, 2 ) 
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s=sin(2*pi*(fs/Iength(a)).*(1:z)); 

data=data1.*s; 

ffx1=m(data(1:64)); 

ffx2=fft(data(65:128)): 

ffx3=fft(data(129:192)); 

ffx4=fft(data(193:256)): 

Z1=ffxr*ffx1; 

Z2=ffx2'*ffx2; 

Z3=ffx3'*ffx3; 

Z4=ffx4'*ffx4; 

ZC=zeros(N,N): 

coherent case Z=Z1 +Z2+Z3+Z4; 

cc=atan2(imag(Z),real(Z)); 

figure,subplot(223),plot(abs(Z)),grid 

title('scatter plot') ylabel('magnitude') xlabel('fl') 

subplot(224),plot(sum(abs(Z))),grid 

titlei'sum plot') xlabel('fl') ylabel('magnitude') 

subplot(221),mesh(abs(Z)),grid 

xlabel('fr) zlabel('magnitude') ylabel('f2') title('mesh') 

subplot(222),contour(abs(Z)),grid 

xlabel('fr) ylabel('f2') title('contour') 

%figure,subplot{311),plot(abs(cc)),grid 

%title('BPSK (coherent case)') 

%subplot(312),contour(cc),grid 

%ylabel('phase') 

%subplot(313),mesh(cc),grid 
%ylabel('phase') 

%xlabel('Z') 

%lncoherent Case 

ZP=abs(Z1)+abs(Z2)+abs(Z3)+abs(Z4): 
figure,subplot(411 ),plot(ZP),grid 
title('BPSK Incoherent (Power) Case') 
subplot(412),plot(sum(ZP)),grid 
subplot(413),plot(sum(ZP')),grid 
ylabel('magnitude') 
subplot(414),contour(ZP),grid 
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